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O ' Abstract 

d : 

We present a general approach to the analysis of gauge stability of 3+1 formulations of 
General Relativity (GR) . Evolution of coordinate perturbations and the corresponding pertur- 
bations of lapse and shift can be described by a system of eight quasi-linear partial differential 
equations. Stability with respect to gauge perturbations depends on a choice of gauge and a 
\ background metric, but it does not depend on a particular form of a 3+1 system if its con- 

^ • strained solutions are equivalent to those of the Einstein equations. Stability of a number of 

known gauges is investigated in the limit of short-wavelength perturbations. All fixed gauges 

CI 

except a synchronous gauge are found to be ill-posed. A maximal slicing gauge and its parabolic 
extension are shown to be ill-posed as well. A necessary condition is derived for well-posedness 
of metric-dependent algebraic gauges. Well-posed metric-dependent gauges are found, however, 
to be generally unstable. Both instability and ill-posedness are associated with perturbations of 
physical accelerations of reference frames. 
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1 Introduction 



Physical analysis of many problems of general relativity (GR) requires a solution of a full 
set of the non-linear Einstein equations. Due to a complexity of these equations, in most 
cases this can be done only numerically. This paper is concerned with gauge stability of 
3+1 approaches to numerical integration of the Einstein equations. 3+1 means here any 
approach in which the equations are split on constraint and evolution parts. Constraints are 
satisfied on an initial three-dimensional space-like hypersurface. Initial data is then evolved 
in time by solving a Cauchy problem. 

The original ADM 3+1 formulation uses a three-dimensional metric 7^ and an extrinsic 
curvature as unknown functions |T[|. In other 3+1 formulations, such as hyperbolic, 
conformal formulations and their combinations, the ADM equations are extended by intro- 
ducing additional variables such as spatial derivatives of 7^, traces of 7^ and K^, conformal 
factors, by forming new combinations of these variables, or by modifying equations with the 
help of constraints [2-18]. These modifications change the nature of the equations so that 
they can become more stable, and the integration can be prolonged [19-23]. Nonetheless, 
calculations (of black hole collisions, in particular) often suffer from instabilities. A stable, 
accurate, long-term integration of the Einstein equations remains as an outstanding problem 
of numerical GR. 

Difficulties with time-integration can be numerical and analytical. Numerical difficulties 
arise when an unstable finite-difference scheme is used to integrate a particular 3+1 set of 
equations. Constructing a stable numerical scheme for a chosen 3+1 set of equations should 
not be a problem if the equations are well-posed and stable. Numerical stability will not 
be discussed in this paper ( see, e.g., |16|, [25[ |26| for numerical aspects). We believe that 
analytical difficulties arising from the nature of the Einstein equations themselves present a 
much more serious and still unsolved problem. 

First of all, a 3+1 set of GR equations must be well-posed [^, p8"] . It is impossible to 



guarantee for an ill-posed system a convergence of solutions when numerical resolution is 
increased. It is known, for example, that a harmonic gauge leads to a symmetric hyperbolic, 
well-posed system [EDI, RUf. Some gauges lead to ill-posed systems (see below). The property 



of well-posedness is local and time-dependent. Using a gauge that maintain well-posedness 
everywhere and at all times is crucial for a long-term stable integration. 

Maintaining well-posedness, however, is not enough. A well-posed non-linear system 
may still have rapidly growing, diverging, or singular solutions so that numerical integration 
will become problematic. Here we call this an analytical instability. Analytical instabili- 
ties can be separated on three types. First, the instability may be related to a physical 
nature of space-time. For example, one can encounter a true singularity where curvature 
invariants become infinite. Second, the instability may be related to violation of constraints 
during time-integration (constraint instability). These may be the energy and momentum 
constraints, as well as additional constraints arising from introduction of additional vari- 
ables into a system. If constraints are satisfied initially, they are automatically satisfied at 
later times. However, a Cauchy problem possesses a much broader class of solutions than 
constraint-satisfying solutions. Consequently, a small initial perturbation may lead to a 
rapid deviation from a constrained solution during integration. The third type is the gauge 
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instability. In a 3+1 approach, a coordinate system is "constructed" during time integration 
according to a pre-determined choice of gauge and initial conditions. Initially small coordi- 
nate perturbations may lead to a divergence of coordinate systems as time progresses. For 
a long-term integration to be successful, a gauge must be both well-posed and stable. 

Little is known about stability properties of 3+1 systems of GR equations in general. The 
analysis is complicated because many different modes are present simultaneously, stability 
is determined by complex non-linear terms in the Einstein equations, and it depends on a 
particular solution and on location in spacetime. Investigation of well-posedness is somewhat 
easier since it requires analyzing only a principal part of the system. Several classes of 
hyperbolic, well-posed 3+1 formulations have been constructed (e.g., |5], [1^, [Tj| |19|), but 



numerical experiments show that these systems are often unstable. There are indications 
that, at least partially, the instability is related to constraint violating modes (e.g., |13|, |I"9|). 



It is long known that a synchronous gauge is prone to the formation of coordinate singu- 
larities (or caustics) . Consequences of this for a numerical integration were discussed 
most recently in |33| . An attempt to characterize the development of coordinate pathologies 
(shocks) in hyperbolic 3+1 formulations has been made in fl3"l| . Perturbations of coordinates 
have been separated, in linear approximation, and studied for a synchronous gauge (see 
Problem 3 in Paragraph 95 of f3T[ and also To our knowledge, instability of gauges 



other than a synchronous gauge have not been analyzed. 

In this paper we develop a general approach to the analysis of gauge instabilities. In 
Section 2 we derive a system of eight quasi-linear partial differential equations that describe 
coordinate perturbations for arbitrary gauges and in arbitrary 3+1 system. We illustrate 
our approach by investigating stability of various gauges in Section 3. A physical meaning 
of gauge instabilities is discussed in Section 4. Our conclusions are given in Section 5. 

2 Equations of gauge perturbations 

Consider a space-time described in a certain coordinate system x a by a four-dimensional 
metric g ab . In what follows we use letters a — h to denote four- dimensional indices 0, 1, 2, 3, 
and letters j-mto denote three-dimensional spatial indices 1, 2, 3; time t = x°. We write 



an interval as 35 



ds 2 = g ab dx a dx b = -(a 2 - ^)dt 2 + 2p i dtdx i + j ij dx i dx j (2.1) 

where a is a lapse, /3, is a shift, (3 l = Y^Pj, Y m lmk = S l k , %j is a metric on a three-dimensional 
space-like hypersurface, and the 4-metric is 

/ -(a 2 + I3if3 l ) flA ab I -^2 ^2 ] be cr . .)■)> 

9ab — a. , 9 ~ 8? iA ° aim , 9ab9 ~ a . {2.2) 



We want to separate, in a linear approximation, the behavior of coordinates from phys- 
ical evolution of a space-time. Under an infinitesimal coordinate (gauge) transformation 
i/j a (x°, x , x 2 , x 3 ) of a four- dimensional coordinate system, 

x a -> x a + V a , (2.3) 
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where ip a are sufficiently smooth infinitesimal functions, the metric transforms as |3l[ 

9ab -> £aft + Sgat , (2.4) 

where 

Sgab = -(V b Va + V a <0&) (2.5) 

are the deviations of the metric, and V denotes a covariant derivative. 

For a metric ( |2.1| ), the deviations Sg a b can also be related to corresponding deviations of 



lapse, shift, and three-dimensional metric, 

U = 6a, Vi = 5f3i, W i3 = 5 li3 (2.6) 

as 

5g 00 = -2aU + 2^(3^ - f k PAW jk 

dgoi = V, (2.7) 
Sgij = Wij 

Substitution of (|2.7|) into (|2.5|) gives us U, Vi and in terms of ip a , 

1 (p&Wq oUr diPo 



v = »[^-PV<--£ + r a i, e ) , (2. 

+ (-) 

where T c ab are four-dimensional Cristoffel symbols. They can be expressed in terms of a, 
Pi, and jij using (|2.2|). Equations ( |2.8| ) and ( |2.9| ) can also be rewritten as a system of four 
quasi-linear partial differential equations for ip a , 

?h. = a u- p% - + {r c 00 + p'p%) ^ c , 

«/ ( 2 - n ) 

a ~ 1 at* + ' 

The deviations ip a , and the corresponding deviations U, V, were arbitrary up to this point. 

During the integration of a 3+1 system in time, certain conditions are imposed on a and 
Pi which specify the choice of gauge. In general, a gauge can be specified by a set of four 
partial differential equations with respect to lapse, shift, three-dimensional metric, and their 
partial derivatives, 

Fa V ' a ' A ' a^' '"' 7iJ ' a^' ") = ' (2 ' 12) 
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where ... in ( j2.12 ) indicate higher order derivatives of a, and 7^. Equations ( 2.12 ) must 



be considered a part of a 3+1 system which is integrated in time along with the rest of the 
equations for metric components, extrinsic curvature, etc. By varying (|2.12|) with respect 



to a, /3i and 7^ and making use of Q2.10 ) we obtain a system of four additional quasi-linear 



partial differential equations relating U and Vi to coordinate perturbations ip a , 
dF a dF a dU dF a d 2 U 

dF a dF a dVi dF a d 2 Vi 

~dfc i + dx 1 + d ( dx b dx c + "' + (2.13) 

dF a f^ c dip A dF a d ( <% 



Together, ( |2.11| ) and ( [2.13] ) form a system of eight equations describing the evolution of 
coordinate perturbations and associated perturbations of lapse and shift in time. Perturba- 
tions of a three-dimensional metric which correspond to evolving coordinate perturbations 
are given by ( |2.10| ). 

At this point, we succeeded in separating the behavior of gauge perturbations from other 
possible perturbations of the solutions of a 3+1 system. Gauge instabilities can be studied 



by investigating a set of eight quasi-linear partial differential equations ( 2.11Q and (|2.13|) 



for coordinate perturbations ip a and associated perturbations of lapse and shift, U and V{. 
During the derivation, we did not use any specific assumptions about a 3+1 system. The only 
assumption implicit to the derivation was the equivalence of constrained solutions of a 3+1 
system to these of the Einstein equations. Therefore, ([2.11 ) and (|2.13| ) must be applicable to 



any such system. In a linear approximation, the behavior of gauge modes of perturbations 
depends only on a choice of gauge and on an unperturbed metric g^. 

Since coefficients of ( |2.11|) , (|2.13j ) are functions of g a b and thus of x a , exact solutions of 



( |2.11|) , (|2.13|) can be found only in some special cases. A further simplification is possible in 
a high-frequency limit where we consider coordinate perturbations on much shorter scales 
compared to those of a base solution g a b. Let us designate 

z(x a ) = {r,U,V i } T (2.14) 

a combined vector of unknowns entering our gauge perturbation equations. We will consider 
in this paper the gauge equations ( |2.13j ) which involve 7^ and its first-order time derivatives 
(cases with higher order derivatives can be treated similarly). Then the resulting gauge 
perturbation equations will contain only first-order time derivatives of z, but, of course, they 
can contain higher-order spatial derivatives of z. We write ( |2.11|) , (|2.13|) symbolically as 

^ _ (0) M , (1) M i dz s , (2) M ij d 2 z s 

where ^Ai, k = 0,1,2,... depend on an unperturbed solution. In the limit of short- 
wavelength perturbations, we look for solutions of (|2.15| ) in the form 

z = C(t)ex V (-lqe k x k ) , (2.16) 
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where I = y— 1, g« = qti is a wavevector, is a unit vector, and g is an absolute value of g^. 

Substituting ( [2.16| ) into (|2.15|) gives an ordinary differential equation for ( (we retain only 
the time-dependence of ^ k 'Ai), 

§ = M(f)C, (2.17) 

with 

M rs (t, a, q) = {0) M rs (t) + ^MUt)leiq + i2) Mi j s (t) ei e jq 2 + ... (2.18) 
When the time-dependence of M can be neglected, the solutions become 

z (x a ) oc exp(u s t - lqe k x k ) , (2.19) 
where wavenumbers u s (q,ei) are determined by the dispersion relation 

det (M - ujN) = , (2.20) 



and N is a unit matrix. If a timescale of the variation of coefficients in M becomes comparable 
to |u;| _1 , the behavior of perturbation with time will not be exponential any more. However, 
a time derivative of z at t = will still be determined by the corresponding value of u. In 
particular, the growth rate of the perturbation amplitudes will be given by Re(u s ). 

An information about Re{uj s ) can be used, first of all, to probe ill-posedness of a gauge. 
If Re(w) — > oo when q — ► oo for at least one uj n and at least in one direction e^, it would 
be possible to construct a harmonic perturbation such that at t = both the perturbation 
and any finite number of its derivatives will be less than any predetermined small number 
e << 1 but will become greater than any large predetermined number A » 1 after a finite 
period of time. This will mean, of course, that the gauge is ill-posed. If, on the other hand, 
the gauge is well-posedQ, than Re(u; s ) will determine the rate of growth of instabilities with 
time. 



3 Gauge stability 

In what follows, it will be convenient to distinguish between the three types of gauges. 
A gauge is called "fixed" if none of F a in ( |2.12| ) involve 7^, and F a (x b , a, Pi) = 0. We 
assume that ( |2.12| ) can be inverted and a, Pi can be expressed explicitly as functions of 



four-coordinates, 

a = a(x a ) } f3i = {3i(x a ). (3.1) 

A fixed gauge does not change when coordinates are perturbed, J^- = Jj^- = 0, and, as a 
result, U — Vi — 0. Equations ( |2.11|) with U — Vi — describe the evolution of coordinate 
perturbations for fixed gauges. 

1 The absence of short-wave harmonic solutions with arbitrary large real increments is a necessary condition for 
well-posedness. Sufficient condition would require proving a continuous dependence of solutions on initial conditions 
for the entire 3+1 system. 
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A gauge is called "algebraic" or "local" if it can be expressed as a function of local values 
of jij and its derivatives, F a (x b , a, Pi, 7^-, ...J = 0, or 

a = a(x a , 7<i ,^,...) , P l = P l (x a ,^,^,...) ■ (3-2) 

Fixed gauges are, of course, a particular case of algebraic gauges. Expressions for U and 
for algebraic gauges reduce to 

u=2 ^( 1-*, _ »*A + 2 * « f iw. - ^ V -. 



V cte a y 



(3.3) 



Finally, a gauge is called "differential" or "non-local" if it cannot be reduced from ( |2.12|) 



to an algebraic form. An algebraic gauge can always be expressed in a differential form by 
simply differentiating (|3.2|). Values of a and Pi for a differential gauge can only be expressed 
as a non-local functional of 7^. 

3.1 Fixed gauges 

To illustrate the approach outlined above, we begin with gauges fl3.1|) which are functions 
of spatial coordinates and time only. A synchronous gauge a = a(t), Pi = is an important 
member of this family. We will show below that it is the only well-posed fixed gauge. All 
other gauges (|3.1| ) are ill posed. 

For a fixed gauge, (|2.11|) become a first-order linear system of partial differential equations 

= + ct. , 

dt dx> /g 4) 

C a = T a 00 + D^ , (3.5) 
D ij = p l p j , (3.6) 



where 



and 



El = 2T a 0k . (3.7) 
For synchronous gauges ( |3.4j ) becomes 

di/j 1 da . 

a (3 8) 
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A general solution of (|3.8| ) is 



ipi 



■ a(t)f , 

7« (f - 



dfo_ 

dx k 



ay 



jh dt 



(3.9) 



where fo{x h ), f l (x k ) are smooth arbitrary function of spatial coordinates (see Problem 3 in 
Paragraph 97 of |3l| ). According to (|3.9|) , there is a continuous dependence of solutions 
on initial conditions in the class of functions with continuous first derivatives. Thus a syn- 
chronous gauge is well posed. 

Next, consider solutions of ( |3.4| ) in a short-wavelength limit. The matrix M in ( |2.18| ) 
becomes 



M = M + IqMt 



(3.10) 



with 



fc° c 1 c 2 c 3 \ 

El E{ El El 

El E\ E\ E\ 

\E\ E\ E\ E\) 



/0 D lk e k D 2k e k D 3k e k \ 



M 1 



ei 
e 2 
\e 3 





















(3.11) 



and ( p.20|) becomes a quartic equation for 



0J, 



(3.13) 



u 4 + d 3 u 3 + d 2 u 2 + d x uj + d = , (3.12) 

with coefficients being polynomial functions of q, 

do = d 0fi + d 0A q + d 0:2 q 2 , 
d\ = d lfi + d 1A q + d 1>2 q 2 
d 2 = d 2 , + d 2 ,i q + d 2)2 q 2 , 

<^3 = ^3,0 , 

where d it j are functions of 7^, a, Pi, and e^. In what follows, we need expressions for 

d 0t2 = (E* m Ei - E^El) D km e k e, + \ - E}E{) d 2 . 2 , 

d x . 2 = ElD^aej - E k k d 2)2 
d 2 ,x = -1(0' 
d 2 . 2 
^3,0 



T\ki\ 



-C° - Et 



and 



(3.14) 



From ( p,12| ) and ( |3.14|) it is clear that the asymptotic behavior of roots of ( |3.12|) with 
q — > 00 must be lo = ^kQ k ^ m , with n, m integer. Values of m, n and coefficients w k can be 
determined by substituting a power series expressions for u into (|3.12| ) and requiring that 
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terms with same powers of q cancel out. The result depends on whether some of dij are zero 
or not. 

If e^2,2 7^ 0, the asymptotic behavior of the roots is 

uj = Uiq + uj + O ( - ) , (3.15) 



with 



Ml,2,3,4 = {0,0, ilV^} , (3.16) 



If o?2,2 = but dip 7^ 0, the asymptotic behavior changes to 

u 2 q 2 ^ + cutq 1 / 3 + u + O (-^J , (3.18) 
where 

Mi.2,3,4 = Sign(d lj2 ) • |rf li2 | 1/3 • {0, 1, ~ ±1^} (3.19) 
If both c? 2 ,2 — ^1,2 — 0, the asymptotic behavior changes again, this time to 

u = u 1 ^ 2 + u + o(4b) , (3-20) 



with 



( ^ = i ,W^±Mr (321) 

The asymptotic behavior is c<j ~ 0(1) if all three coefficients c?2,2 = ^1,2 = ^2,1 = 0. In 
particular, for a synchronous gauge a = 1, = 0, the increments are determined by 
det(~f km ^ - u8f) = 0. 

For fixed algebraic gauges with zero shift we have D %: > = 0, c?2,2 — ^1,2 — 0, and from 
( |3.21|) the asymptotic behavior is 



^1,2= ±g 1/2 (l-lT^)\/^ + 0(l), «m = 0(1). (3.22) 



Cj in (|3.5| ) reduces for $ = to 



= n = -a 7 ^ . (3.23) 
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Since is a unit but otherwise arbitrary vector, there always will be harmonic solutions with 
Re(u) ~ g 1 / 2 and the gauge will be ill-posed unless lapse is spatially constant, a = a(t). 

Consider now gauges with fa 7^ 0. According to ( |3.16|) , Re(w) ~ 0(1) for all wavevectors 
not orthogonal to shift, faei 7^ 0, so that these perturbations will not cause ill-posedness. 
However, it can be seen from (|3.17 ) that for small fa « 1 the increment Re(uj) ~ ±/3 _1 C l ej 



can become arbitrary large. Thus, a mode of perturbation unstable for fa — cannot be 
eliminated by applying a small shift in the direction of its propagation. 

For wavevectors which are orthogonal to shift, faei = 0, we find — {E\D k i — 
E^D % i)eiej = 0, and the asymptotic behavior is again given by ( |3.22| ) but now with C l 
determined by the general formula ( |3.5|) . If C l fa 7^ 0, all orthogonal harmonic solutions will 
have Re(uj) ~ g 1 / 2 , and the gauge will be ill-posed. We now show that this is always the 
case, that is, C % and fa cannot be co-linear. Expressing T% c in terms of three-dimensional 
quantities, we find that C % can be written as 

C = 7*(^ + +fl*)ff. (3-24) 

Thus, co-linearity of C % and fa requires that fa must satisfy a system of partial differential 
equations 

dfa d(a 2 - 7 mn p n p m ) _ 

~dt dx~ k " ' [ } 

These equations depend on 7™ m and, thus, fa must depend on ■j nm as well. This contradicts, 
however, to our initial assumption (Oj) that a and fa are explicit functions of spatial coor- 
dinates and do not change when r f im are perturbed. It is easy to see that the only solution 
of ( |3.25p independent of r f %n is fa = 0, and that this solution is possible if a = a(t). We 
finally conclude that, unless shift is zero and lapse is spatially-constant, there always will be 
initially arbitrary small harmonic solutions of (3.11) which can be made arbitrary large after 
a finite period of time by an appropriate choice of That is, among fixed algebraic gauges 
( |3.1|) only synchronous gauges fa = 0, a = a(t) are well-posed. All other gauges (|3.1|) are 
ill-posedQ. 

3.2 Algebraic gauges 

As a next example, consider gauges with metric-dependent lapse and fixed shift, 

a = a{x a , llJ ) , fa = fa{x a ) (3.26) 

From ( p.3|) we obtain 

U-tp-fa*,-^) , V i = 0, (3.27) 



dlu V 13 dxi 



2 In a one-dimonsional case where both a metric and its perturbation are dependent on one spatial coordinate 
x = x , we have e = {1, 0, 0}, and a non-zero shift along the ^-coordinate leads to (5id 7^ and to a well posed gauge. 
This however, is not a contradiction since coordinate perturbations orthogonal to gauge were not allowed. 
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and the coordinate perturbation equations (|2.11| ) become 



dt 

dijji 



da 2 \ dtp 



d^fij J dx 



+ T a 00 + fa^ + 



da 2 



d% 



'j 



■pa 

ij 



(3.28) 



dt dx' 1 

Comparing gauge perturbation equations Q3.28| ) and ( |3.4D we observe that (|3.28|) becomes 



identical to 



if the definition of D tJ is changed from (3.6) to 

da 2 



D 



ij 



faP J + 



d'ji 



(3.29) 



a 2 

In what follows, we assume that Jp- and , thus, D lJ are symmetric. We obtain a dispersion 

^^^^ 

relation similar to ( |3.12| ), and, after taking a limit of q — > oo, get an asymptotic behavior of 
wavenumbers u which dependence on coefficients dij is similar to that described by equations 
( |3.15| - ^2T1) . Now, however, D 1 ^ ^ fa 1 ft and thus g? 2 , 2 = D^e^tj may be both positive, zero, 
and negative. For g? 2j2 < we have two roots with Re(uj) ~ q, and the gauge is ill-posed. If 
d 2j2 = 0, d% t 2 7^ 0, we see from ( |3.19| ) that one root will have Re(u) ~ g 2//3 and the gauge 
is again ill-posed. If d 2)2 — 0, d\^ = 0, d 2 ,i 7^ 0, then one root has Re(u>) ~ q 1 ! 2 and the 
gauge is ill-posed as well. Having g? 2 2 = rfi )2 = rf 2 ,i = and, thus, the asymptotic behavior 
uj ~ 0(1) for all four roots is impossible. A reasoning similar to that used for fixed algebraic 
gauges above shows that this requires fa to be functions of the metric and contradicts to 
the assumption ( |3.26|) . We finally conclude that gauges (|3.26|) are ill-posed unless D^e^ej is 
strictly positive, i.e., unless 



dji 



eiej > 



(3.30) 



is satisfied for every e^. 

One simple gauge of type ( |3.26| ) which has been used in a variety of works is a gauge 
which depends on the determinant of a three-metric, 



a 



a(x a ,~f) , 7 = det(jij) 



For this gauge 



D 



ij 



cry 



ij 



(dry 



yfidjij), and, since 7 y is positive definite, the gauge will be ill-posed unless 

da 2 



$7 



> . 



(3.31) 



(3.32) 



(3.33) 



More general metric-dependent gauges may be investigated in a similar way. For example, 
gauges with both a and fa functions of the metric, 



a = a(x a , 7^) , fa = fa(x a , 7y) 



(3.34) 
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have both U =fi and V t ^ 0, 



07u 



(3.35) 



and the resulting coordinate perturbation equations become more complicated 



dipp 

dt 

dipj 

dt 



dx i 



dxi 



dx k 



(3.36) 



Here 



D' 



Olij Olij 



^k — Zi Ok ~ Zi i] 



(3.37) 



9li 



The matrix M remains the same, but M\ now becomes 



Mi 



/0 L> lfc e fc D 2k e k D 3k e k \ 

ei G^e* Gfe* Gf e fe 

e 2 G X 2^k G^ek G^ek 

Ve 3 G* fc e fc Gfe* G| fe e fey / 



(3.38) 



The dispersion relation is still a quartic equation ( |3.12| ) but its coefficients now are forth- 
order polynomials in q. The analysis of the roots ^1,2,3,4 is a more complicated algebraic 
problem. In general, ( |3.34| ) may be either ill- or well-posed depending on the functional form 
of a and 

It must be stressed here that a well-posedness is not a guarantee of gauge stability. Metric- 
dependent gauges satisfying (|3.30|) are generally unstable with the increment of instability 
Re(u) = uq given by formulas ( p. 17] ) with coefficients determined according to (|3.5| ), (|3.7|) , 
( |3.14|) , and (|3.29|) . In particular, it can be shown using these formulas that the gauge (|3.31|) 
can be unstable for a wide range of background solutions. When ( |3.34| ) is well-posed, it can 
be both stable and unstable depending on a particular background solution. 

A class of gauges often considered in the literature is || 



dt H dx l Jy ' v • 



(3.39) 
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where f(a) is an arbitrary function. It was found that / > is a necessary condition for 
hyperbolicity of first-order 3+1 formulations of GR introduced in that paper. The gauge 
( |3.39|) is, in fact, equivalent to an algebraic gauge || 

Vl = F(a) (3.40) 

with / = aF 1 , and contains as its members such gauges as a harmonic slicing (/ = 1) 
and a "1+log" slicing (/ = I /a). One can see that the condition / > derived from the 
analysis of hyperbolicity of the entire 3+1 system of equations is equivalent to our condition 
of well-posedness W- > derived from the analysis of gauge modes alone. 



A more general family of first-order 3+1 systems has been derived in |L9[ using a gauge 



log(a<T ff ) = Q{x a ) , Pi = fo{x a ) . (3.41) 

This gauge belongs to a family of algebraic gauges (|3.26| ) as well. It was found in 
that having a metric-dependent, densitized lapse with a > is a necessary condition for 
a hyperbolicity of 3+1 systems considered in that work. Again, it is easy to see that the 
requirement cr > is equivalent to the condition of well-posedness derived in this paper from 
the analysis of gauge instabilities. 

3.3 Differential gauges 

As an example of a differential gauge, we consider a parabolic extension of a well known 
maximal slicing gauge Y^ViVja = K^K^a , $ = 0, |37H , 

f)a 1 

_ = - (-yVViVjOt - K^K ija ) , A = , (3.42) 

where e > is a constant. We follow the same general procedure as that used above for fixed 
and algebraic gauges. Expanding covariant derivatives in Q3.42 ) and taking into account 



that, for Pi = 0, the extrinsic curvature = —i^Ttij, we rewrite (|3.42|) as 



dt ~^ ^ dx % dxi ^ ij dx k Act > & > 

and find the derivatives of F , 

dF 7r kl 7r kl dF kl ■ dF _ l3 8F 

"7 Ah , oO \ 7 



9° ' 4 « 2 ' 9(B) " ~ " ' d(oSh) 1 ' dm " ' 

nr A — — 7 7 &w + < 7 7 A « + 7 7 A «) ft? • (3 - 44) 

^ = "o?7W T +7 7 -7 7 j , and 



di^fs) 2dx nK ' ' 11 ' ' ' ' 07^ " 2a 
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Substituting (|3.44 ) into (|2.13|) , and combining (|2.13|) with ( |2.11|) we obtain the system of 
five quasi-linear partial differential equations 



where 
C k 



all + r^ c , 
dip 



dt 

dipj 

dt 
dU 



dx i 



+ 



Aa 2 



a dx l dxi 



2B i j>k f^lH . + cl ^Po + Dlj d^ 



dxWx k 



dxi 



dxj 



71 



—it + 2— T° 0i + 25^4 , D ij 
a J a J 



ik 



214^ 



a 



ir ij /dT c - 



9 ^ 



2A ij + 2B lk > j T\ k , and 

dT c - 
+ 2B ij ' k —^ . 
ox k 



The dispersion relation becomes 



(3.45) 



(3.46) 



det 



/ 




r 1 

1 00 


r 2 

1 00 


r 3 

1 00 


« 




2^ + Ie l9 


21^-0; 


2F 2 

Zi 01 


2r 3 

Zi 01 







2r° 2 + Ie 2 g 


Zi 02 


2r 2 2 - u 


2r 3 

Zi 02 







2r° 3 + Ie 3 g 


Zi 03 


2V 2 

Z,L 03 


2r 3 3 -^ 





V 


M 40 


M 41 


M 42 


M 43 


M 44 - 



(3.47) 



where 



M 40 

M 44 : 

M 44 



-\J- lC l e iq 



TV ^ GiGj 2 

q 

a 



- (J 1 - \D i] e jq + 2B l3 ' l ej eiq 2 ) , and 



The dispersion relation ( J3.47| ) is a fifth-order equation with respect to u, 

—to 5 + d^to 4 + d 3 uj 3 + d 2 uj 2 + diuj + d = , 



(3.48) 



(3.49) 
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where 



g?4 — d^fi + d/^ ig + d^ 2 q 2 , 

d 3 = 4,o + d 3A q + d^ 2 q 2 , 

d 2 = d 2 ,o + d 2 ,\q + d 2)2 q 2 + d 2 , 3 q 3 

di = di t0 + d 1A q + d 1<2 q 2 + <ii )3 g 3 

do = d ,o + d 0jl q + d 0j2 q 2 + d 0t3 q 3 



(3.50) 



and coefficients dij are functions of an unperturbed metric g a \,. In what follows, we need 
explicit expressions for 



"4,2 

d^i 
d 2 ,3 



I (2e ^fl^e^e* + t^Ime*) = -Irf 4 , 2 ro e fc 



t ^ft i nk 



(3.51) 



From ( [3.49 ), Q3.50|) we obtain the following asymptotic behavior of the roots. 
U j 1 = -d 4 2 q 2 - d A)1 q + 0(1) , 




+ 1 =± 1-1 



Q 



1/2 



roo e fcl 



+ 0(1) 



(3.52) 



Since 7 U is positive definite, we have d^ 2 > 0, and thus Re(ui) < in the limit of q — > oo. 
This root does not cause an instability. However, one of the roots u 2t3 has a positive Re(u) 
growing with q as 



Reiui) oc g 



1/2 



roo e fc 



9 



1/2 



da 



dx k 



Since the direction of 6, is arbitrary, a maximal slicing gauge will be ill posed unless a is 
spatially constant. The result we obtain for a maximal slicing and its parabolic extension is 
exactly the same as that we found earlier for fixed algebraic gauges (see ( |3.22|) ). It is easy to 
see that in both cases the asymptotic g 1//2 -behavior of the roots comes from the presence of 
r oo V?j term in the right-hand side of the first equation in either ( |3.4| ) or ( 3.45 ). Differentiating 
first equation with respect to time and substituting second equation into it, we obtain in 
both parabolic-like equation for ipo 



d 2 ^ 



J- r 



+ other terms 



8t 2 ' 00 dx i 

which is ill-posed if initial conditions are specified at t — const, (see also next section). 
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4 Physical meaning of gauge instabilities 



4.1 Gauge transformations in a flat spacetime 

To illustrate our analysis of gauge stability, we must consider now the entire evolution part 
of the Einstein equations and compare its stability with that of the system ( |2.11| ), ( |2.13| ). 



We will use an ADM 3+1 system with metric 7^ and extrinsic curvature Kij as unknown 
variables, 



dig 
dt 



3 ' 



= -ViVja + a fa + KK %] - 2K im Kf) 
+ p m V m K tJ + K im V 3 f3 m + K jm Vi(r , 

and will consider a simple class of one- dimensional solutions of ([4.1] ) 

7^ = diag(y, 1, 1) , ^ = diag(K, 0, 0) (4.2) 

with 7 and dependent on t and x = x 1 only. In vector notations, (|4.1|) becomes 

(4.3) 



<9t dx 



where 



u=(^), A=l_^_^ 9 J^ M , fl= Kas^aKg g) • (4-4) 



One can easily verify that all components of the Riemann tensor are identically zero for 
any metric (|4.2|) satisfying ( |4.3| ). That is, ( |4-2p describes a flat spacetime in non-Galilean 
coordinates^]. One can also verify that the constraint equations are automatically satisfied 
for any metric ( f4.2|) . Therefore, ( f4.3|) can be unstable only with respect to gauge instabilities. 
Eigenvectors and eigenvalues of B in (|4.3|) are 



A = ± ^u=( 4 ' 3 «-°- ) ) , (?) • (4.5) 



3 The only nontrivial component of the Riemann tensor in this case is i?ioio- Since choice of a, (3 is a choice of 
a coordinate system, it is sufficient to show that JJioio = in a coordinate system with a = 1, f3 = 0. In this case, 
non-zero Cristoffels are 

F i _ r i _ 1 £7 r o _ 1 37 pi _ 1 a 7 
1(11 %~dt ' in "29t' lll_ 2^' 

and 

Rowi - -2~w~ + 9llT ^ -~2W + ^{m, 
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For (3 7^ 0, there is a complete set of of eigenvalues and eigenvectors, so that (|4.3|) is strongly 
hyperbolic and well posed. For /3 = 0, a complete set does not exists, and the system reduces 
to a parabolic equation 

<9 2 7 { ada\ d'y 9 K 2 da d 2 a 

' - )-^ + 2a 2 2— K + 2a— - (4.6) 



dt 2 \ 7 dx J dx 7 dt dx 2 

da _ d 2 a 
dx dy 2 



which is ill-posed (a parabolic equation |? = jc4 is ill posed if initial conditions are defined 



at x — const, e.g., |p8fl ). 

Analysis of the principal part is not sufficient, however, for investigation of insta- 
bilities of (|4.3|) . We need to carry out a full stability analysis of this system. A linearized 
version of ( f4.3|) is 

- = C(u)u + B(u)- (4.7) 
where u is now an unperturbed base solution, u is a vector of perturbations, and 



d _dA dBdu _ £g -2a 

du du dx ~ \ aR2 - — 1 — - 2/3 — - + 4E_/3 ^ 25^ _ P 97 _ 2aJC 

y 7 2 cte 27 2 9x t 2 " <9x 7 2 9z 7 3 cte 7 9x 7^ 9ie 7 

(4.8) 

For harmonic solutions 

u oc exp(c<jt — Igx) , (4.9) 
in the limit q — > oo we then obtain a dispersion relation 

\\C -IqB -Iuj\\ =lu 2 + d 1 u + d = , (4.10) 



with 



and 



7 jdxj ' (4.11) 
do = do,o + Mo,i9 + rf o,2<? 2 , 



(3 d 7 \ /2<9/? _ 2(3 3-1 _ 2aK\ 
7 2 cte / V 7 <9x -f 2 dx 7 / 



2a (^K 2 da 1 57 /3 2i*Td/3 4^/5^7 

\ 7 2 <9x 27 2 <9x 7 2 dx 7 2 <9x 7 3 <9x 
2(3 ( 18(3 (3 d-f aK\ a da 

"0,1 = o o 7T" + 



(4.12) 



do,2 

T 



7 dx 'j 2 dx 7 / 7 9a; 
/? 2 ^ 
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For (3 7^ 0, rf 0) 2 7^ and we get an asymptotic behavior 

w = ±lq ^. ± ^^ ± l^l_^£ + i?l T (^ + l?E) +0 (l ) (413) 

7 2(3 dx 7 2 <9x 7 7 <9x \ 7 7 dx ) \q ) 

We see that although the system is well-posed, it has Re{ui) ~ 0(1), and is unstable with 
respect to small perturbations (unless If = = |f = 0, and K > 0). If |^ 7^ 0, the 
increment of the instability Re(tu) ~ (3 — > 00 when /? — ^ 0. 

For (3 = and I s 7^ 0, we have d 0:2 = 0, do,i 7^ 0, and the asymptotic behavior is 



« = ±(1 + 1)^1 + 0(1) (4.14) 

which shows that Re(u) ~ g 1 / 2 and the system is ill-posed. The conclusions agree with the 
results of gauge stability analysis of the previous section. 

4.2 Instability of a synchronous coordinate system 

Synchronous coordinate systems are formed by acceleration-free test particles. It is well 
known that in such systems the metric determinant \ will vanish after a finite time because 
the time-lines of a reference frame will necessarily intersect one another on a certain caustic 
hypersurface f3~IJ . This is correct in both flat and curved spacetimes. On a caustic, one of 



the principal values of the metric tensor vanishes, whereas the corresponding contravariant 
component tends to infinity. In the vicinity of a caustic, using an arbitrariness in the selection 



of spatial coordinates one can write a general four-metric as [32 



goo = -l, g i = 0, g np = ln P = a np , g n3 = 7 nl = r a nl , g n = 711 = r a n (4.15) 

where r = t—x l , indices n,p take values 2, 3, and are non-singular functions of x 2 , x 3 , and 
t. Coefficients are connected by a single relation which is a consequence of the Einstein 
equations (see Appendix in [0). Equation ( (4.15Q determines a general "quadratic" charac- 
ter of approaching a caustic in a synchronous gauge. In ( |4.15| ), some of the contravariant 
components of g ab tend to infinity as r -2 . Because of that, the process described by (|4.15|) 
is sometimes called a "blow up" instability. 

It is important to distinguish between the formation of a caustic (blow up) and the 
instability of a synchronous gauge with respect to small perturbations considered in this 
paper. To clarify this important point, let us consider the following example. Consider 
one-dimensional solutions of ( [4.3| ) for a = 1 and (3 = 0. The equations ( |4.3| ) become 



dt ' dt 7 ' v ; 

and a general solution of ( |4.16|) isQ 

7 = (A(x) + B(x)t) 2 , K = -B(x)(A(x) + B(x)t). (4.17) 

4 Physical meaning of A and B in this solution is discussed in the Appendix 
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If at some point we have A > 0, B < 0, a caustic will form and both 7 and K will become 
zero at some future time t c = —A/B > 0. Note that the process of caustic formation may 



be global if = or local otherwise. Using (|4.17|) it is easy to write the deviation of 7, 

SIS db function of perturbations of A and B, 

57 = 2(A + Bt)(A5A + BSBt) , (4.18) 

where 5 A = AA/A, SB = AB/B, and 8 A, 5B << 1. This is a weak instability described by 
a power function of t instead of an exponential function^]. Both the formation of a caustic 
and the instability are due to the presence of initial velocities of the reference frame with 
respect to an inertial Galilean frame. However, the caustic is due to the intersection of 
converging trajectories of test particles, whereas the instability refers to a divergence of test 
particles from unperturbed time-lines. 

4.3 Gauge instabilities in accelerating systems 

Let us now discuss a physical meaning of ill-posedness and instability of non-synchronous 
gauges. Consider first a Rindler reference frame |38| where this meaning is most clear. Let 



T,X l ,X 2 ,X 3 be a Cartesian coordinate system in a Minkowski space-time with an interval 

ds 2 = -dT 2 + (dX 1 ) 2 + (dX 2 ) 2 + (dX 3 ) 2 (4.19) 

The Rindler frame is constructed from the world lines of particles moving with accelerations 
along an X-axis of this coordinate system in such a way that the accelerations are constant 
in proper time of the frame. Such a frame is "rigid" in the sense that It does not deform 
in proper time. For an accelerating frame to be rigid, accelerations of the particles must 
be different with respect to the inertial frame. If accelerations of particles with different 
X 1 and same time T are equal, the distance between particles measured in the inertial 
coordinate frame will stay constant. Due to a Lorentz contraction, distances measured in 
the accelerating frame will depend on the velocity of the frame, and will change with time. 

The Rindler frame can be described using coordinates x, t (we omit two other spatial 
coordinates for brevity) as 

T = xsinh{gt), X 1 = xcosh{gt) , (4.20) 

where g is a constant. In these coordinates the interval becomes 

ds 2 = -{gxfdt 2 + dx 2 , (4.21) 

and it is evident that the geometry of the co-moving space is this frame is time-independent. 
Coordinate lines of t, T 2 = X 2 — x 2 , are world lines of uniformly accelerated particles. An 
acceleration a measured by a co-moving observer in an orthogonal coordinate system x a is 



3S 



a 



\lli^a k , a* = P 00 / 5oo • (4.22) 



A formal stability analysis gives in this case a growt h rat e of perturbations u) — —K/j comparable to that of an 
unperturbed solution. As it was discussed in Section 3, (p.l9|) is not valid in this case. 
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For a system with acceleration along the X 1 -axis we obtain 

1 da 



and for a Rindler frame 



dx 



(4.23) 



a R = - (4.24) 

x 

In a Rindler frame, the coordinate x is equal to a physical distance I = J* y/^fndx from the 
point x = where physical acceleration is infinite. Therefore, ( |4.24j ) can be written in terms 
of quantities which have a direct physical meaning, 

a R = y . (4.25) 

The dependence ( |4.25|) of a physical acceleration on a physical distance provides the rigidity 
of a Rindler frame. 

The perturbation equations ( |4.7| ) for a Rindler frame can be reduced to a single equation 

where 7 is a perturbation of 711. Consider the following exact solutions of ( |4.26| ) (modes of 
perturbation) : 

1. 7 = const « 1. This corresponds to a uniform change of the length scale along the 
x-axis. Such a perturbation changes both the acceleration a = -(1 + 7) -1 ^ 2 and the 
physical distance I = J^(l +j) 1 ' 2 dx = (1 + ^j) l ^ 2 x in such a way that still a = l/l. The 
condition ( [4.24 ) is not violated and the frame remains rigid. 



2. 7 — At, At << 1. The perturbation does not depend on a? and describes a deformation 
of the frame with constant velocities with respect to a Rindler frame. Formally this is 
an instability but of the type similar to that existing in synchronous gauges. As in the 
previous case, it can be shown that the relation a = l/l holds. 

3. 7 = x~ a2 ^ 2 (ci sinh(v4t) + C2 cosh^t)). This time the initial perturbation depends on x 
and grows exponentially. Suppose c 2 > so that 7 > 0, then |2 < 0. The acceleration 
is again a(x) = ^(l + 7(a;))~ 1 / 2 . The physical distance, however, is I = L (l+*f) l ' 2 dx = 
x(l + 'j(x')) 1 ^ 2 with some x' < x. As a result, a(x) > l//(x) is now greater than that 
required to maintain rigidity, and the deformation will grow with time. A perturbation 
7 = x A 1 9 (ci sin(y4t) + c 2 cos(At)) is also a solution of ( |4.26| ), but now |2 < and this 
particular perturbation is stable. 

We see that the reason for the instability is a mismatch between the deformation of the 
frame and the prescribed acceleration of particles which leads to further frame deformation. 
We can imagine, for example, a deformation with 7 = everywhere for x less than some 



20 



value xo, 7 = 70 = const > at x > x\ > xq and 7 increasing linearly between xq and x\. 
The difference between accelerations at xq and X\ will be finite, but by tending x\ — ► so we 
can make 7(^/7(0) between the neighboring points so and xi to increase with any desired 
rate. This corresponds to the mathematical property of ill-posedness of fixed gauges with 
spatially non-uniform lapse. It is also obvious that when a and 7 depend on t and x, the 
instability will mean a deviation from an already deforming reference frame. 

In one-dimensional case, any metric with non-zero shift can be obtained from a metric 
with (3 = by a transformation of a time coordinate t = i(t, x). Such a transformation does 
not change the accelerations of test particles which remain "chronometric invariants" 4^ . 



However, due to a non-zero shift, perturbations are now "advected" , that is, their coordinates 
change continuously whereas accelerations are still prescribed at fixed coordinates. As a 
result, j(x) cannot grow with an arbitrary rate. An apparent growth of the instability at 
fixed x will decrease with increasing f3. In our stability analysis, this corresponds to a changes 
from an ill-posed to a well-posed but unstable gauge. 

4.4 Instability of rotating reference frames 

In one-dimensional cases considered above the cause of instability and ill-posedness was the 
perturbation of accelerations of test particles forming the reference frame due to perturbation 
of positions of these particles. In more than one dimensions, there is another physical factor, 
rotation and Coriolis forces associated with it. However, as we will see below, this factor does 
not change the nature of instability. Consider, for example, a uniformly rotating reference 
frame 

a = 1, pi = {-ttx 2 , p 2 = ttx 1 , 0}, 7^ = diag(l, 1, 1) . (4.27) 
where Q is an angular velocity. For (|4.27|) we have 



& = r 00 + pp k r Jk = r 00 = {oV, n V, 0} , (4.28) 

and, according to our gauge stability analysis ( |3.23| ), a perturbation with a wave vector 
Qi — Q^i will grow with the increment 

Re{w) = Qq 1/2 vVei + x 2 e 2 . (4.29) 

Note, that C % are proportional to components of physical acceleration eij in ( |4.23| ). The 
increment oc q 1 ^ 2 indicates that the gauge is ill-posed. The reason is the same as in a one- 
dimensional case. It is related to a radial acceleration and not to Coriolis force. Note that 
at the axis of rotation x 1 = x 2 = we have C l = and there is no instability (Re(u) = in 
( 4.29|) ), but of course, there is a Coriolis force. Thus we conclude that in a three-dimensional 



case with rotation the main physical reason for gauge instability and ill-posedness is the 
same acceleration a 1 . 

5 Conclusions 

In this paper we presented a general approach to the analysis of gauge stability of 3+1 for- 
mulations of GR. Gauge modes of perturbations can be separated, in a linear approximation, 
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from other modes of perturbations and studied independently. A system of eight quasi-linear 
partial differential equations (|2.11|) , (|2.13| ) describes the evolution of coordinate perturba- 
tions ip a and the corresponding perturbations of lapse and shift, U and with time. The 
gauge stability depend on the choice of gauge and on an unperturbed four-metric g a b, but it 
does not depend on a particular form of a 3+1 system of GR equations. 

Well-posedness and stability of several gauges was investigated. We demonstrated that 
all fixed gauges, i.e., gauges that are functions of coordinates only, are ill-posed with the 
exception of a synchronous gauge a — ot(t), /3j = 0. This gauge is well-posed, but it is prone 
to the formation of coordinate singularities (caustics) and it is unstable. 

It is known that maximal slicing tr(Kij) = prevents the formation of coordinate singu- 
larities by applying accelerations to test particles forming a reference frame in such a way 
that |2 cx tr(Kij) = and the local volume element remain 7 1 / 2 = const. [[HJ]. This, and a 



singularity-avoiding properties of maximal slicing allow in many cases a much longer numer- 
ical integration than using a synchronous gauge. However, both the maximal slicing gauge 
and its parabolic extension ( |3.42| ) are ill-posed due to the presence of acceleration-related 
unstable modes. Computations using these gauges will blow after a long enough period of 
integration. 

Stability of metric-dependent algebraic gauges has been investigated as well. In partic- 
ular, the necessary condition of well-posedness of gauges with metric-dependent lapse and 
fixed shift ( |3.26|) was formulated. In addition to the formation of caustics, algebraic well- 
posed gauges with spatially dependent lapse and shift are susceptible to instabilities caused 
by perturbations of accelerations in deforming reference frames. The reason for the unstable 
behavior of these hyperbolic gauges is the same as that for ill-posedness of fixed, parabolic, 
or elliptic gauges. However, due to an "advective" property of hyperbolic gauges, the growth 
rate of acceleration-related unstable modes becomes wavelength-independent, Re{u) ~ 0(1), 
in the limit q — > 00. 

An investigation of stability of 3+1 formulations of GR in this paper is limited to an 
investigation of gauge instabilities. By studying ill-posedness and stability of (|2.11|) , ( |2.13|) 
one can tell if a 3+1 sets of GR equations using a particular gauge will be ill-posed or 
unstable. All gauges studied in this paper were found either ill-posed or unstable. None of 
the gauges investigated in this paper are suitable for a long-term stable integration of GR 
equations. Gauges with better stability properties must be found. 

Stability of a gauge does not mean, however, that a 3+1 system using this gauge will be 
stable. As was mentioned in the introduction, another source of ill-posedness and instability 
is associated with violation of constraints. Constraint instabilities do depend on a particular 
form of a 3+1 system. The analysis of hyperbolicity in jn| provides an illustration of 
this statement. A part of the hyperbolicity conditions derived in this work (their equation 
(2.36)) does not involve gauge at all, and is dependent on how the constraint equations are 
incorporated into the system. 
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Appendix: Physical meaning of A and B in ( 4.171) 
Let us introduce the following: new coordinates t, 

T _ t + Ux 1 



xl g + Ui [ x fdU\ ydy (A.1) 



X 2 = x\ X 3 = x 3 , 

where T, X 1 are Cartesian coordinates in a Minkowski spacetime g a b = diag(—l, 1, 1, 1). In 
these new coordinates we get the metric with the components 

5-00 = -1 , 9u = {A + Bt) 2 , g 22 = #33 = 1 , others g ik = , (A. 2) 

where 

It is clear that (A.l) is a generalization of a Lorentz transformation with the velocity U 
depending on x 1 . 
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